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ABSTRACT 

When  moisture  and/or  temperature  are  suddenly  changed  on  the  boundary  of  a 
solid,  stresses  and  strains  are  introduced  and  they  can  be  further  aggravated  by 
the  presence  of  stress  raisers  such  as  voids  or  cavities.  A  time  dependent  fi¬ 
nite  element  procedure  is  developed  for  solving  the  hygrothermal  stresses  around 
a  circular  cavity  in  a  finite  plate.  Numerical  results  are  displayed  graphically 
for  the  T300/5208  graphite  fiber-reinforced  epoxy  resin  material.  The  size  of 
the  hole  relative  to  the  plate  is  varied  for  three  different  cases  such  that  the 
interaction  of  moisture  and  temperature  is  investigated  in  conjunction  with  changes 
in  the  solid  geometry. 

Possible  failure  sites  are  also  examined  by  application  of  the  strain  energy 
density  criterion.  These  locations  are  determined  from  the  stationary  values  of 
the  strain  energy  factor.  The  hygrothermal  influence  tends  to  move  the  failure 
site  away  from  the  cavity  while  the  mechanical  load  gives  the  opposite  effect. 

The  proportion  of  the  energy  stored  by  hygrothermal  and  mechanical  disturbances 

is  investigated. 
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INTRODUCTION 


The  general  influence  of  moisture  and/or  temperature  on  the  stresses  and 
displacements  in  laminated  composite  materials  [1,2]  has  received  increased  at¬ 
tention  in  recent  times.  It  is  known  that  when  a  laminate  absorbs  moisture, 
its  mechanical  stiffness  and  strength  are  degraded  and  recovery  is  incomplete 
after  desorption.  Heat  can  also  degrade  a  material.  These  environmental  influ¬ 
ences  can  interact  so  that  the  stress  state  of  the  material  is  dependent  on  both 
temperature  and  moisture  in  its  surroundings.  A  theory  of  diffusion  which  in¬ 
corporates  the  interaction  between  temperature  and  moisture  can  be  found  in  [3]. 
Phenomenological  arguments  leading  to  coupled  equations  governing  the  simulta¬ 
neous  diffusion  of  moisture  and  heat  were  further  elaborated  in  [4].  All  the 
physical  models  led  to  the  same  system  of  differential  equations  although  the 
coefficients  related  to  the  basic  thermodynamic  properties  of  the  solid  differed. 
Discussed  were  the  reciprocal  effects  of  heat  and  moisture. 

The  stresses  produced  in  a  plate  by  the  hygrothermal  strains  associated  with 
the  diffusion  processes  described  earlier  have  been  calculated  [5].  Both  the 
conditions  of  suddenly  applied  temperature  and/or  moisture  on  the  plate  surface 
were  considered.  Since  the  temperature  and  moisture  concentration  in  the  plate 
vary  with  time,  the  stresses  also  fluctuate  and  tend  to  zero  when  the  temperature 
and  moisture  concentration  become  uniform.  The  situation  when  the  moisture  dif¬ 
fusion  coefficient  is  temperature  dependent  was  treated  in  [6]  for  the  symmetric 
boundary  conditions  which  produce  no  bending.  Coupling  of  moisture  and  heat  was 
found  to  be  inherent  in  the  case  of  transient  temperature  boundary  condition  for 
a  given  moisture  content.  Depending  on  the  magnitude  of  the  surface  temperature 
change,  the  stresses  predicted  from  the  coupled  and  uncoupled  theory  can  differ 


anywhere  from  20  to  80%.  Bending  is  produced  when  the  boundary  condition  is 
skew-symmetric  [7]. 


The  present  investigation  is  concerned  with  the  moisture,  temperature,  and 
stress  fields  for  the  problem  of  a  plane  body  containing  a  circular  hole.  The 
coupled  diffusion  equations  with  polar  symmetry  are  first  solved  by  a  time  de¬ 
pendent  finite  element  procedure.  Since  the  elastic  deformation  is  assumed  not 
to  be  coupled  with  moisture  and  temperature,  the  stresses  can  be  solved  indepen¬ 
dently  once  the  diffusion  field  is  determined.  Two  types  of  transient  boundary 
conditions  are  considered.  They  are  the  sudden  change  of  temperature  and  mois¬ 
ture  on  the  circular  hole.  Numerical  results  are  displayed  graphically  and  dis¬ 
cussed  in  connection  with  the  minimum  strain  energy  density  criterion  for  loca¬ 
ting  possible  failure  sites.  The  accuracy  of  the  time  dependent  finite  element 
procedure  developed  in  this  work  is  tested  by  solving  the  one-dimensional  diffu¬ 
sion  problem  of  a  plate  whose  solution  is  known  in  closed  form  [5],  The  results 
agree  very  well  and  are  given  in  the  Appendix. 

FINITE  ELEMENT  FORMULATION 

A  time  dependent  two-dimensional  finite  element  method  will  be  developed  to 
solve  the  coupled  diffusion  equations  of  heat  and  moisture: 

Dv2C  -  (C-XT)  =  0 

(1) 

Pv2T  -  (T-vC)  =  0 

in  which  v2  =  32/3x2  +  a2/ay2  stands  for  the  Laplacian  operator  in  two  dimensions 
and  t  is  the  time.  In  equation  (1),  T  is  temperature  and  C  is  the  mass  of  mois- 


i  » 
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ture  per  unit  volume  of  void  space  in  the  solid.  The  diffusion  coefficients  D 


and  V  have  units  of  area  per  unit  time,  and  the  coupling  coefficients  \  and  v 
have  units  of  mass  per  unit  volume  per  unit  temperature  and  the  reciprocal,  re¬ 
spectively.  These  equations  are  relatively  easy  to  solve  in  the  one-dimensional 
case  and  when  the  coefficients  are  constant  and  boundary  values  of  temperature 
and  moisture  content  are  held  constant  between  occasional  moments  of  sudden 
changes  [5]. 

Referring  to  Figure  1,  the  boundary  value  problem  to  be  considered  here  is 
inherently  two-dimensional  and  multiply-connected  with  an  inner  boundary  and 
outer  boundary  r~.  The  enclosed  region  is  denoted  by  R.  For  t<0,  the  tempera¬ 
ture  and  moisture  fields  are  such  that 


T(x,y,t)  =  TQ(x,y) 
C(x,y,t)  =  CQ(x,y) 
while  for  t>0  they  change  to 


(2) 


T(x,y,t)  =  TQ(x,y)  +  AT(x,y,t) 
C(x,y,t)  =  CQ(x,y)  +  AC(x,y,t) 


(3) 


Since  equations  (1)  cannot  be  solved  generally  by  analytical  means,  a  finite 
element  procedure  will  be  developed. 

Bcuic  ^onmuZcutlon.  In  order  to  apply  a  scheme  used  in  variational  calculus,  the 
following  scalar  functions  <j>i  and  <j> ^  are  introduced: 


ii  k* CM 
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(4) 


*1  ■  V  'I  £<f>2  +  <f>2J  +  c  It  (C-'T))dxdy 

*2  ’  V  '!  [<|7>2  +  (f)2’  +  T  It  <T-vC)ldxdy 


The  desired  solution  to  equations  (1)  can  be  obtained  by  requiring 

<5$1  =  54>2  =  0  (5) 

Let  the  body  in  Figure  1  be  divided  into  m  triangular  elements  with  n  nodes.  The 
moisture  C  and  temperature  T  at  the  nodes  will  be  denoted  with  the  subscripts 
i,j  and  k  while  C  and  T  will  refer  to  the  values  in  the  element  A.  For  a  linear 
relation,  C  and  T  can  be  written  as 


(i  *  h  c<Yk  -  Vj>  *  (J,j-yk)x  *  <vxj)yJ 


with  a  being  the  area  of  a  typical  element  shown  in  Figure  1.  Note  that  the  ex¬ 
pressions  for  Nj  and  may  be  obtained  by  the  cyclic  permutation  of  subscripts. 
A  system  of  linear  equations  can  thus  be  obtained  with  the  help  of  equations  (5) 
and  (6): 
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3<j>, 

3C7  =  0,  i  =  1 ,2,—  ,  n 

3<t>o 

YC  =  0,  i  ■  1,2,—  ,  n 


Let  the  vector  C  and  T  be  defined  as 

**  0+ 


rcii 

Tr 

C2 

T2 

. 

• 

• 

II 

•  [ 

• 

• 

C 

T 

_n^ 

_n_ 

(9) 


Applying  the  conditions  in  equations  (8)  to  (4)  yields  the  expressions 


HC  +  DKC  -  XHT  =  0 

HT  +  PKT  -  vHC  =  0 

^  »  *>  'X  A»  •• 


(10) 


where  dot  represents  differentiation  with  respect  to  time.  The  elements  of  the 

★ 

matrices  H  and  K  for  the  finite  element  A  are  given  by 


★ 

If  the  origin  of  the  coordinate  system  (x,y)  in  Figure  1  is  placed  at  the  cen- 

(A) 

troid  of  the  element  A,  then  H:-'  becomes 

'  vl 


H‘j’  ■  is  [Vj  *  T?  b1bj<x1+xj+xk>  +  T2  (bidj  *  bjdi><Vl  +  Yj  +  Vk> 


*  Ti  didj(yi+yj+yk)] 


-6- 


i 


V 


HSj}  =  fi?  y  (ai +  bix  +  v)(aj +  bjx  +  v)dxdy 

Ki{})  -  h  (bibj  +  W 


oi) 


in  which 


-  xjyk 


-  x 


kyj 


bi  " 


yj'yk’ 


di  " 


'  X; 


02) 


and  the  remaining  quantities  in  equations  (11)  can  obviously  be  obtained  by  the 
cyclic  permutation  of  the  subscripts.  The  quantities  H . .  and  K, •  can  be  obtained 

J  ^  J 

by  summing  up  H^,  and  K^,  for  all  elements  A,  B,  etc., 

i.e. , 


H.  .  =  +  h!!^  + 

ij  ij  ij 


Kij 


K(  A) 
Kij 


k(b) 

Kij 


(13) 


Vzcompoiing  matriazi  and  vector.  In  order  to  solve  equations  (10)  numerically, 
it  is  convenient  to  rearrange  it  to  the  forms 


(|  -  v)HC  +  y  KC  +  VK1  »  0 
(-  -  A)HT  +  -  KT  +  DKC  =  0 

'  ****  yj  A***  ***** 


(14) 


The  time  portion  of  the  problem  will  be  solved  analytically  so  as  to  give  an  ade¬ 
quate  treatment  of  the  transient  nature  of  the  boundary  conditions. 
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Equations  (14)  will  be  re-structured  by  decomposing  the  matrices  and  vectors 


into  parts  referring  to  the  inner  region  with  subscript  I  and  to  the  boundary 
nodal  points  with  subscript  8.  This  leads  to 


(t  -  v)HtCt  +  Y  KtCt  +  PKtTt  =  -  4  -  v}HTpCrj  -  T  *<toCd  -  PKtcTd 
X  X  ~1  'X  ~IB~B  X  ~IB~B  ^IB^B 


(|  -  X)HtTt  +  -  KtTt  +  DKtCt  =  -  (1  -  X)HtoTr  -  -  KtdT0  -  DKtdCr 
'X  v  v  ~IB*B  v  ^IB^B 


(15) 


As  t  becomes  infinitely  large,  the  steady  state  condition  is  recovered  and  equa¬ 
tion  (15)  reduces  to 


KTC  T 


+  JSibSb 


JSiIi 


+  -ibIb  = 


(16) 


In  what  follows,  the  transient  boundary  conditions  of  sudden  moisture  and/or  tem 
perature  change  will  be  considered. 

BOUNDARY  CONDITIONS 

Let  CQ  and  7^  be  the  initial  values  of  the  moisture  and  temperature  vectors 
while  Cf  and  T^  are  those  corresponding  to  the  final  values.  The  quantities  ACg 
and  ATg  represent  the  increment  change  of  these  vectors  on  the  surfaces  Tj  and 
r j j  as  shown  in  Figure  1,  i.e.. 


Cr  =  C  +  ACn 


If  lo  +  aIB 


(17) 
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Sadden  MotiiuAe  Change.  The  first  type  of  transient  boundary  conditions  in¬ 


volves  a  sudden  change  of  the  moisture  condition  on  Tj  while  the  temperature  re¬ 
mains  constant.  No  changes  occur  on  j .  This  can  be  expressed  as 


AC 


L  0  , 


on  r 

on  r 


I 

II 


AT  =  0,  on  r,  and 


(18) 


Under  these  considerations,  equation  (15)  becomes 


(1  -  '■IBiti  * 
<7  -  * 

Taking  the  Laplace 


x  J£i£i  +  %Ii  ■  •  l  ~iB~f " 

?  ill,  +  OKlC,  -  -  DKIBCf  - 
transform  of  equations  (19) 


"SibIo  -  (T  ■  v)BiB4£fs(t> 
?  SibIo 

such  that 


(19) 


II 


£l 


it  can 


=  /  TT(t)  exp(-st)dt 
o 


(20) 

oo 

=  /  C j ( t)  exp(-st)dt 
o 

be  shown  that  by  eliminating  Cj,  the  following  expression  is  obtained: 


■'i 

I 

■i 


(s2M  +  SN  +  A)(Tj-T0)  =  -  (y  -  v)(HIKI1  +  HIBACf) 
in  which  M,  N  and  A  are  given  by 

^  AT 
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(21) 


i 


M  *  1  (1  -  v)(l  -  XJH.K^H. 

~  D  X  v  -I 


(22) 


N  *  TT  “  1)Hr 
~  D  Xv  ~I 


Equation  (21)  may  be  solved  as  an  eigenvalue  problem  with  a  solution  of  the  form 


Tt  -  T  =  l  e.f.(t) 

~l  -o  -^i  ~i  1 
J  1 


(23) 


where  n  is  the  number  of  nodes.  The  eigenvectors  e^  obey  the  relation 


(A  -  ufM)e1  =  0,  i  =  1,2,...,  n 


(24) 


For  a  nontrivial  solution,  the  determinant  |A-u?M[  must  vanish  which  may  be  evalu¬ 
ated  to  yield  a>- .  The  scalar  function  f..(t)  is 


f  (t)  ,  y*P(-W> 

^  /y^-1 


sinh(a)it/vf-l ),  !  y1  I  >  1 


(25) 


In 


equation  (25),  and  y.  stand  for 


“1  *  -  <1  -  '■>*!  <HiS'Kib4£b  *  Sib4£b> 


(26) 


y  -  =  e^Ne  . 

1  2o).  ~1  ~  -1 


Similarly,  the  same  procedure  can  be  repeated  to  eliminate  Tt  in  the  expres- 
sions  for  the  Laplace  transforms  of  equations  (19).  This  leads  to 
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(27) 


-11- 


(32) 


(v  ■  mil  +  v  ISiXj  +  °Kj£i  *  -  y  ISielf  '  d-ib^o  ’  C  '  x^iBA'B5^t^ 

Following  the  procedure  discussed  earlier  where  use  was  made  of  Laplace  trans¬ 
form,  it  is  found  that 


II  -  lo  ■  4, 


(33) 


The  function  h . ( t )  is  found  to  be 


★ 

6, 


h i ( t )  =  [1  -  exp(-Y..it)cosh(.it/^  + 


exp(-Yio)it) 


“i  ’'M'"’ 


*  M 


* 

(<i  -  --■)  sinh(u.  t/yX-T)  for  |y..|  >  1 


(34) 


in  which  8i  and  are  given  by 


•;  ■  -  i7>  sisre 


*1  *  D  *’  '  L* 


(35) 


Returning  to  equations  (32),  Cj  -  ^  can  also  be  evaluated  and  it  takes  the  form 


CT  -  C  =  l  e<k.(t) 

~i  ~0  -j =]  "“1  1 


(36) 


with  k^(t)  being  given  as 


a.exp(-y.ui.t)  _ 

k,.(t)  *  - - - sinh((i)^t/y^-l )  for  |y,.|  >  1 


i  Ti 


(37) 
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such  that 


«*  -  -  §  -  ‘'sIsa'^Ib  (38> 

This  completes  the  formation  of  the  time-dependent  portion  of  the  problem. 

The  geometric  portion  will  be  solved  numerically  by  the  method  of  finite  element 
in  the  numerical  examples  to  be  followed. 

NUMERICAL  EXAMPLES:  DIFFUSION  PROCESS 

The  geometry  to  be  treated  is  that  of  a  square  plate  LxL  containing  a  circu¬ 
lar  hole  of  radius  a.  Referring  to  Figure  2,  if  the  boundary  conditions  on  the 
hole  or  Tj  is  independent  of  9,  then  the  problem  possesses  1/8  -  symmetry.  Only 

the  shaded  area  needs  to  be  analyzed.  Considered  are  three  different  L/a  ratios 

★ 

with  L  being  equal  to  8,  16  and  1,352  .  This  corresponds  to  progressively  weaker 
interaction  of  the  hole  with  the  plate  boundary  as  the  hole  radius  a  is  always 
kept  at  unity.  The  distance  b  =  j  -  a  is  kept  constant  so  that  the  numerical  re¬ 
sults  for  all  three  cases  can  be  compared  on  the  same  graph.  The  finite  element 
grid  patterns  are  shown  in  Figures  3  to  5  and  they  are  self-explanatory. 

The  plate  material  is  made  of  an  epoxy  resin  used  for  the  T300/5208  graphite 
fiber-reinforced  composite.  The  coupling  constants  were  determined  in  [6]  and 
they  are  D/O  =0.1,  x  =  0.5  and  v  =  0.5.  Numerical  results  for  the  moisture  and 
temperature  distribution  around  the  circular  cavity  are  obtained  for  the  condi¬ 
tions  specified  by  equations  (18)  and  (31)  which  will  be  discussed  separately. 

Mo-ostuAe  Change..  For  Case  I  with  L=8,  a  plot  of  the  normalized  moisture  change 

(C-Co)/(Cf-CQ)  versus  (r-a)/b  is  given  in  Figure  6  for  different  values  of  the 

★ 

In  this  case,  L  represents  the  diameter  of  a  circular  plate. 
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time  parameter  Dt/b2.  For  small  time  t,  the  moisture  concentration  drops  rapidly 
as  a  function  of  the  distance  (r-a)/b  and  only  the  material  near  the  circular  hole 
boundary  is  affected.  As  time  increases,  the  material  away  from  the  hole  is  also 
influenced  by  moisture  change  and  the  decrease  in  (C-CQ)/( Cf-CQ)  becomes  more 
gradual.  The  variations  of  temperature  with  distance  measured  from  the  hole  are 
exhibited  in  Figure  7.  The  temperature  peaks  near  the  hole  for  small  time  and 
their  values  tend  to  decrease  and  move  away  from  the  hole  as  t  is  increased.  The 
results  for  Case  II  with  L=16  are  not  appreciably  different  from  those  shown  in 
Figures  6  and  7  and  hence  will  not  be  displayed.  When  L=1  ,352  such  that  the  hole 
diameter  is  decreased  considerably  in  size  as  compared  with  that  of  the  plate,  a 
significant  change  of  results  are  observed.  They  are  referred  to  as  Case  III 
illustrated  in  Figures  8  and  9.  The  moisture  changes  in  Figure  8  for  different 
times  are  seen  to  take  place  only  in  the  material  close  to  the  hole  boundary. 

The  temperatures  in  Figure  9  do  not  peak  as  significantly  as  those  in  Figure  7 
when  the  hole  is  closer  to  the  plate  edge.  A  comparison  of  results  for  all  the 
three  cases  is  made  in  Figures  10  and  11  in  terms  of  the  average  moisture  and 
temperature  defined  as 

ca»e  '  V  /  “*•  Tave  ’  7  /  ™  <39> 

Displayed  in  Figure  10  is  ( Cave-CQ )/ (C^r-C^ )  versus  »^t/b.  As  it  is  to  be  ex¬ 
pected,  a  significant  uptake  in  moisture  as  a  function  of  time  is  observed  for 
Case  I  where  the  circular  hole  occupies  a  greater  portion  of  the  plate.  The 
presence  of  the  hole  becomes  less  and  less  significant  as  its  size  is  reduced 
going  from  Case  I  to  Case  III.  Similarly,  the  average  temperature  peaks  at 
fOtl b  =  0.55  for  all  cases  with  the  largest  peak  corresponding  to  Case  I  when 
the  hole  and  plate  edge  interaction  is  the  strongest.  This  is  shown  in  Figure 
11. 
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TempestatuAe  Change..  If  the  temperature  on  the  hole  boundary  Tj  is  suddenly 
changed  from  T-  to  while  the  surface  moisture  is  kept  constant.  Figure  12  dis¬ 
plays  the  drop  in  temperature  as  a  function  of  (r-a)/b  for  Pt/b2  =  0.011,  0.044 
and  11.111.  The  material  near  the  hole  experiences  more  severe  temperature  drop 
for  small  time  t.  As  more  of  the  material  is  subjected  to  temperature  changes, 
this  influence  becomes  more  gradual.  The  corresponding  changes  in  moisture  C-Ci 
are  similar  to  those  shown  for  T-T..  in  Figure  7  when  moisture  boundary  condition 
is  specified.  In  fact,  the  results  for  C-C^  can  be  obtained  from  those  in  Figure 
7  for  T-T..  by  the  multiplication  factor  P/D  since  x  =  v  =  0.5  for  this  particular 
example.  Therefore,  C-C^  will  peak  near  the  hole  and  then  decrease.  Figure  13 
gives  the  decrease  in  temperature  versus  distance  for  Case  III.  The  trend  is 
the  same  as  that  in  Figure  12  except  that  the  influence  is  confined  closer  to 
the  hole  boundary.  The  values  of  C-Ci  for  Case  III  again  differ  from  those  in 
Figure  9  for  T-Ti  by  the  factor  P/D  and  hence  will  not  be  repeated. 

TRANSIENT  HYGR0THERMAL  STRESSES 

Once  the  moisture  and  temperature  distribution  in  R  are  known,  the  stresses 
can  be  obtained  in  a  straightforward  manner  as  follows: 


a  =  E(r£o)  (40) 

in  which  the  stress  tensor  a  has  the  components  and  a  „  and  the  strain  ten- 

r  9  r9 

sor  has  the  components  er,  Eg  and  erg  referred  to  the  polar  coordinates  r,  ?  and 
z  in  Figure  2.  For  plane  strain,  the  matrix  E  for  an  isotropic  and  homogeneous 
material  is  given  by 
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E(l-0 


i,  yo-y. 


s =  n+vp)o-2y  V°‘V’  u 


0,  (l-2vp)/2(l-vp) 


where  E  is  the  Young's  modulus  and  vp  the  Poisson's  ratio.  Assuming  that  ez  =  0, 
the  transverse  normal  stress  component  az  can  be  found  as 


°2  *  yar+<y  ■  +  S(C-C0)1 


The  coefficient  of  thermal  expansion  is  a  and  of  moisture  expansion  is  e.  The 
hygrothermal  thermal  strain  eQ  in  equation  (40)  takes  the  form 


£o  =  (HV  9o 


in  which  e  is  defined  by 
o 


9o  *  “<T-To>  *  s<c'c0> 


There  remains  the  determination  of  e  owing  to  the  nonuniform  distribution  of  C 
and  T  throughout  the  elastic  plate  containing  a  circular  cavity.  This  will  be 
accomplished  by  application  of  the  finite  element  procedure. 

Let  q  be  the  equivalent  nodal  force  which  is  statically  equivalent  to  the 

A * 

tractions  applied  on  the  element.  It  can  be  expressed  as 


q  =  Qa  -  p 


J  CA-a'j' 


such  that 


Q  =  /  BTEBdV,  p  =  /  BTEe  dV 
~  y  ~  —  y  -  o 

and  a  is  the  equivalent  nodal  displacement  corresponding  to  q: 


The  matrix  B  is  given  by 
B  -  [B  - ,  B.,  Bj 

~  -vK 

in  which 


(46) 


(47) 


(48) 


Bi 


yrV  0 
°*  Vxj 


yj"ym 


(49) 


and  Bj  and  8^  may  be  written  down  by  the  cyclic  permutation  of  indices.  Since 
both  Q  and  p  in  equations  (46)  are  known  quantities,  a  may  be  found  from  equation 

V  *4 

.  ★ 

(45)  and  hence  the  strain  e  is  determined  since 

At 


Note  that  the  displacement  ij  is  related  to  a  as 


u  =  CN. I,  HU,  N,I] 

~  I  ~  J 


f!i 

l*i 


with  I  being  the  identity  matrix  and  N, ,  N-,  etc.,  are  given  by  equation  (7)  and 

~  i  j 


u  has  the  components 

V 


u  =  Vi  *  Vj  *  Vk 
v  =  Vi  ♦  Vj  +  Nkvk 


r 


1 


C  =  [8.,  B.,  B.] 

<■»  1  ~J 


*i 

Sk 


(50) 


It  is  now  obvious  that  equations  (43)  and  (50)  may  be  inserted  into  equation  (40) 
to  yield  the  hygrothermal  stresses. 


For  the  T300/5208  resin  material,  the  following  material  properties  will  be 
used  for  the  hygrothermal  stress  calculations: 


a  =  4.5  x  10‘5  m/m°C 
6  =  2.68  x  10'3  m/m/%  H20 
E  =  3.45  GN/m2  (5  x  105  psi) 
vp  »  0.34 


(51) 


Again,  the  discussion  on  stresses  for  the  cases  of  moisture  change  and  temperature 
change  will  be  presented  separately. 

M o-u-tuAe  Change.  In  all  of  the  cases,  the  stresses  will  be  expressed  in  MN/m2 

and  plotted  against  the  dimensionless  distance  (r-a)/b  such  that  r=a  refers  to 
the  points  on  the  hole  boundary  and  r  =  a+b  =  1/2  refers  to  points  on  the  plate 
edge.  The  dimensionless  time  parameter  Vt/b2  is  varied  from  0.011  to  11.111. 

The  stresses  ar,  a0  and  a2  for  Case  I  are  given  in  Figures  14  to  16.  The  compo¬ 
nent  ar  in  Figure  14  is  zero  at  r=a  and  r  =  1/2  as  required  by  the  free  stress 
boundary  conditions.  It  is  compressive  in  the  interior  of  the  plate  and  varies 
with  time.  The  peak  of  the  compressive  stresses  tend  to  move  away  from  the  hole 
as  time  is  increased.  The  circumferential  stress  component  c9  plotted  in  Figure 
15  is  compressive  near  the  hole  and  becomes  tensile  at  a  finite  distance  r  which 
increases  with  time.  A  similar  trend  is  observed  in  Figure  16  for  the  transverse 
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normal  stress  except  that  the  effect  is  not  as  pronounced.  When  the  hole  di¬ 
ameter  is  small  in  comparison  with  the  plate  width  L,  it  is  interesting  to  note 
from  Figures  17  to  19  referred  to  as  Case  III  that  the  stress  variations  tend 
to  be  confined  closer  to  the  hole  boundary.  The  compressive  portion  of  or  and 
a0  are  greater  in  magnitude  than  those  for  Case  I.  However,  the  tensile  portion 
of  o 0  and  az  greatly  reduced.  Refer  to  Figures  18  and  19  for  Case  III  and  the 
results  in  Figures  15  and  16  for  Case  I.  Since  stresses  for  Case  II  do  not  dif¬ 
fer  significantly  from  those  for  Case  I,  they  are  not  presented. 

TempcAatuAe  C/iange.  If  the  hole  is  subjected  to  a  sudden  change  in  the  surface 
temperature  as  specified  in  equation  (31),  the  stresses  acquire  an  oscillatory 
character  changing  from  tension  to  compression.  The  variation  depends  on  the 
elapsed  time.  Figure  20  shows  that  ar  is  tensile  for  small  time  near  the  hole 
and  becomes  compressive  as  r  increases.  The  maximum  value  of  in  tension  oc¬ 
curs  at  intermediate  time  as  it  becomes  entirely  compressive  for  large  time. 

The  variations  of  a.  and  a,  in  the  material  ahead  of  the  hole  at  different  times 
are  illustrated  in  Figures  21  and  22.  A  plot  of  ag  versus  time  for  r=a  is  dis¬ 
played  in  Figure  23.  It  clearly  shows  that  oa  reaches  a  peak  at  Pt/b2  =  0.127 
and  then  decreases  and  becomes  compressive.  The  component  az  will  have  a  similar 
behavior.  Figures  24  to  26  give  the  stress  results  for  Case  III.  As  the  hole 
size  is  reduced,  the  compressive  portion  of  of  in  Figure  24  tends  to  dominate 
while  the  tensile  portion  is  greatly  diminished.  The  magnitude  of  both  aQ  and 
a2  are  reduced  appreciably  in  Case  III  and  the  results  are  given  in  Figures  25 
and  26.  In  general,  the  elevation  of  the  stress  state  decreases  with  the  ratio 

of  L/a.  Figure  27  gives  a  summary  of  the  values  of  ct0  for  all  three  Cases  I, 

II  and  III.  The  solid  curves  correspond  to  Pt/b 2  =  0.011  and  the  dotted  curves 
to  Pt/b2  =  0.444. 
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FAILURE  CRITERION:  STRAIN  ENERGY  DENSITY  THEORY 


Having  obtained  the  hygrothermal  stresses  a„,  a.  and  a,  around  the  circular 
cavity,  it  is  natural  to  inquire  into  possible  sites  of  failure.  A  criterion 
that  has  been  used  successfully  for  predicting  failure  of  solids  due  to  yield¬ 
ing  and/or  fracture  is  the  strain  energy  density  theory  [9,10].  The  theory  as¬ 
sumes  failure  to  occur  when  the  energy  stored  within  a  unit  volume  of  material 
reaches  a  critical  value.  This  energy  density  can  be  computed  from  the  stresses 
as  follows: 


37  =  [a2  +  +  02 


1+v. 


+  °Z) 


+ 


(52) 


The  location  of  failure  corresponds  to  dW/dV  being  a  minimum  and  the  physical 
meaning  of  this  condition  can  be  best  interpreted  by  resolving  dW/dV  into  the 
sum  of  two  components.  The  first  component  (dW/dV )v  is  associated  with  volume 
change  and  the  second  (dW/dV)d  with  shape  change.  The  locations  of  dW/dV  minimum 
corresponds  to  failure  by  volume  change  and  are  most  likely  to  result  in  fracture 
while  dW/dV  maximum  corresponds  to  failure  by  yielding.  These  are  relative  mini¬ 
mum  and  maximum  values  of  dW/dV  and  occur  exclusively  within  the  material,  not 
including  any  physical  boundary.  They  are  most  conveniently  obtained  by  taking 
derivatives  of  dW/dV  with  respect  to  the  position  angle  e  of  the  radial  vector  r 
measured  from  a  reference  point  to  a  possible  failure  site.  Refer  to  Figure  2. 
Hence,  the  relation  dW/dV  =  S Jr  is  often  used.  The  quantity  S  is  extracted  from 
dW/dV  as  the  1/r  coefficient  and  is  known  as  the  strain  energy  density  factor. 
Numerical  results  of  dW/dV  will  only  be  given  for  Case  I  where  L =8  as  indicated 
in  Figure  3. 
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MoIMuac  Change.  Making  use  of  equation  (52),  the  strain  energy  density 
dW/dV  is  computed  as  a  function  of  the  radial  distance  r  for  different  time  t. 

When  the  hole  boundary  is  subjected  to  a  sudden  change  in  moisture,  the  minimum 
value  of  dW/dV  or  (dW/dV)ml-n  tends  to  increase  with  t  reaching  a  limit  as  time 
becomes  increasingly  larger.  This  is  shown  in  Figure  28.  Failure  is  assumed 
to  occur  when  (dW/dV)^  reaches  the  critical  value  of  dW/dV  or  (dW/dV)c  for  a 
given  material.  The  maximum  (dW/dV)  ..  is  approximately  10.5  x  10~4  MJ/m3.  The 
trend  of  the  curve  in  Figure  28  implies  that  damage  due  to  moisture  boundary  con¬ 
dition  is  a  long-time  effect. 

TempeA&tu/LC  Change.  The  variations  of  (dW/dV )m^n  with  time  exhibit  a  different 
character  when  the  hole  experiences  a  sudden  temperature  change.  Figure  29  shows 
that  there  are  two  sets  of  (dW/dV)min.  One  has  a  larger  peak  (dW/dV)™?*  =  1.75 
x  10“4  MJ/m3  with  Vt/b2  =  0.23  and  occurs  at  a  larger  distance  away  from  the 
hole,  r/a  =  3.15.  The  other  has  a  lower  peak  ( dW/dV )mi n  =  0-75  x  10~4  MJ/m3  with 
Vt/b2  =  0.20  and  occurs  at  a  smaller  distance  away  from  the  hole  r/a  =  1.10.  For 
the  same  time  t,  (dW/dV )min  is  seen  to  be  larger  at  distances  further  away  from 
the  hole.  The  strain  energy  density  criterion  seems  to  suggest  that  failure  due 
to  hygrothermal  stresses  alone  is  more  likely  to  occur  at  approximately  one  di¬ 
ameter  distance  away  from  the  hole. 

Sup^xpoi-ition  ofi  Mechanical  SViebid .  Since  in  most  applications  mechanical  loads 
are  also  present,  it  would  be  natural  to  inquire  into  the  combined  influence  of 
mechanical  and  hygrothermal  stresses.  In  the  case  of  a  circular  hole  of  radius 
a  subjected  to  uniaxial  applied  stress  a  ,  the  stress  field  is  given  by 
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Cr-£[(1  -£)  Ml  COS26] 

ae  =T-U]  +W]  -  (1  cos2e] 

«re  -  -  F-  ^  +  &  '  si"29 


(53) 


Consider  the  supposition  of  an  applied  tensile  mechanical  stress  of  aQ  =  0.2 
MN/m2,  the  top  curve  in  Figure  29  can  be  combined  with  the  (dW/dV)ml-n  obtained 
from  equations  (53).  This  leads  to  the  results  given  in  Figure  30  for  0=0° 
and  90°  which  correspond  to  planes  parallel  and  normal  to  the  direction  of  ap¬ 
plied  stress.  The  value  of  (dW/dV)™9*  is  seen  to  occur  on  a  plane  normal  to  the 
applied  tension  at  a  distance  approximately  r  =  3a.  The  situation  is  reversed 
when  the  applied  stress  aQ  is  compressive.  Figure  31  shows  that  the  most  likely 
failure  site  is  now  in  a  plane  parallel  to  the  applied  load,  i.e.,  9=0°  along 
which  (dW/dV)min  is  larger.  In  general,  the  failure  is  assumed  to  occur  when 
the  first  { dW/dV )mi- n  reaches  (dW/dV)c- 

As  the  magnitude  of  the  applied  mechanical  stress  aQ  is  gradually  increased 
to  1.0  MN/m2,  the  predicted  failure  site  tends  to  move  in  closer  to  the  hole 
boundary  and  the  lower  curve  in  Figure  29  becomes  more  dominant.  For  aQ  =  2.0 
MN/m2,  it  is  seen  from  Figure  32  that  the  predicted  failure  site  is  much  closer 
to  the  hole.  (dW/dV)j^*  for  °0  Positive  and  aQ  negative  both  occur  at  approxi¬ 
mately  r/a  equal  to  1.1. 

CONCLUDING  REMARKS 


The  hygrothermal  stresses  induced  by  the  sudden  change  of  moisture  and/or 

temperature  at  the  boundary  of  a  circular  cavity  are  determined.  A  time  dependent 
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finite  element  procedure  is  developed  in  which  the  time  portion  of  the  problem 
was  solved  analytically  by  means  of  Laplace  transform.  All  calculations  are 
carried  out  on  the  CDC  6400  computer  and  the  numerical  results  are  believed  to  be 
accurate.  An  estimate  of  this  accuracy  can  be  evidenced  by  the  numerical  computa¬ 
tion  of 


_  1  D+V  /I 

Yi  2w.j  D  Xv 


-  1 


e  H.e.  =  2.0083 
"I 


(54) 


which  agrees  with  the  exact  value  obtained  from  equation  (30).  Another  check  can 
be  seen  from  the  solutions  of  the  one-dimensional  slab  problem  which  has  been 
solved  analytically.  This  is  given  in  the  Appendix.  Special  care  has  also  been 
given  to  scaling  the  grid  patterns  for  Cases  I,  II  and  III  where  relative  dimen¬ 
sions  of  the  hole  and  plate  are  varied. 

In  addition  to  determining  the  coupling  effects  between  moisture  and  tempera¬ 
ture,  the  strain  energy  density  criterion  was  applied  to  investigate  possible 
failure  sites.  Several  interesting  results  were  observed.  First,  the  energy 
state  due  to  the  hygrothermal  influence  alone  tends  to  dominate  in  a  region  ap¬ 
proximately  one  diameter  away  from  the  circular  hole  while  mechanical  loading  ex¬ 
erts  more  effects  on  the  stress  and  energy  states  close  to  the  hole  boundary. 
Thus,  the  precise  location  of  failure  will  depend  on  the  combined  influence  of 
hygrothermal  and  mechanical  stresses.  It  should  be  noted  that  the  present  analy¬ 
sis  did  not  consider  coupling  between  diffusion  and  mechanical  deformation.  Such 
an  interaction  will  be  left  for  future  investigation. 
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APPENDIX:  MOISTURE  AND  TEMPERATURE  DISTRIBUTION 


This  section  considers  the  special  problem  of  a  slab  of  thickness  h  which 
coincides  with  the  z-direction.  The  external  surfaces  of  the  slab  being  parallel 
to  the  xy-plane  are  subjected  to  sudden  change  in  moisture  and/or  temperature. 

The  problem  is  one-dimensional  in  space  as  variations  in  C  and  T  occurs  only  as 
a  function  of  z.  The  solutions  for  the  coupled  moisture  and  temperature  diffu¬ 
sion  problem  are  given  to  illustrate  that  the  present  time-dependent  finite  ele¬ 
ment  procedure  yields  the  same  results  as  those  obtained  analytically  [5]. 

The  grid  pattern  is  given  in  Figure  33  and  the  same  constants  D/D  =  0.1,  • 

=  0.5  and  v  =  0.5  are  used  for  the  numerical  computation.  Figures  34  and  35  qive 
plots  of  (C-Co)/(C^-Co)  and  ( T-Tq)/ v( C^-CQ)  versus  2z/h  for  the  case  when  the 
moisture  on  the  slab  surfaces  are  suddenly  raised  from  CQ  to  while  the  surface 
temperatures  are  held  constant.  Both  the  moisture  and  temperature  levels  tend 
to  increase  with  the  parameter  4Dt/h2.  Similar  plots  are  displayed  in  Figures  36 
and  37  for  the  sudden  application  of  uniform  temperature  to  the  slab  surfaces. 

As  mentioned  earlier,  these  results  when  compared  with  the  closed  form  solutions 
show  that  the  finite  element  method  developed  here  is  indeed  reliable. 
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Figure  4  -  Grid  pattern  for  Case  II  (L  3  16  units) 
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Figure  6  -  Normalized  moisture  change  versus  radial 


distance  for  Case  I  with  AC  =  constant 
and  AT  =  0  on  Tj 
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Figure  7  -  Temperature  change  versus  radial  distance 
for  Case  I  with  AC  =  constant  and  AT  =  0 
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Figure  8  -  Normalized  moisture  change  versus  radial 
distance  for  Case  III  with  aC  =  constant 
and  aT  =  0  on  Tj 
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Figure  9  -  Temperature  change  versus  radial  distance 
for  Case  III  with  AC  =  constant  and  AT  = 


Case  I  ( Ls  8) 


Figure  12  -  Variations  of  normalized  temperature  with  radial 
distance  for  Case  I  with  AT  =  constant  and 
AC  =  0  on  Tj 


Figure  13  -  Variations  of  normalized  temperature  with 
radial  di'stance  for  Case  III  with 
&T  *  constant  and  aC  =  0  on  j 
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Figure  15  -  Circumferential  stress  versus  distance  for 
Case  I  with  sudden  moisture  change  on  r. 
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Figure  17  -  Radial  stress  versus  distance  for 

Case  III  with  sudden  moisture  change  on  rT 
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Casem  (  L  =  1,352 ) 
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Figure  18  -  Circumferential  stress  versus  distance  for 
Case  III  with  sudden  moisture  change  on  rT 
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Case  I  (  L  =  8) 

AT  =  constant-,  A C  =  0  on  Tj 

Figure  20  -  Radial 
Case  I 

stress  versus  distance  for 
with  sudden  temperature  change 
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Figure  23  -  Circumferential  stress  as  a  function  of  time  for 
Case  I  with  sudden  temperature  change 


Figure  24  -  Radial  stress  versus  distance  for 
Case  III  with  sudden  temperature 
change  on  rT 
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Figure  26  -  Transverse  normal  stress  versus  distance 
for  Case  III  with  sudden  temperature 
change  on  rT 
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Figure  28  -  Minimum  strain  energy  density  as  a  function  of  time  for 
Case  I  with  aC  =  constant  on  rT 
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Figure  30  -  Variations  of  (dW/dV)^  with  time  for  Case  I 
AT  =  constant  and  applied  tensile  load 


-56- 


cr0  ( Compression) 


0.01  0.10  1.00 

2>t/b2 

Figure  31  -  Variations  of  (dW/dV)  -n  with  time  for  Case  I 
AT  =  constant  and  applied  compressive  load 
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Figure  36  -  Variations  of  temperature  with 
distance  for  slab  subjected  to 
&T  =  constant  and  AC  =  0 
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Figure  37  -  Variations  of  moisture  with 
distance  for  slab  subjected 
to  aT  =  constant  and  AC  =  0 
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Figure  31  -  Variations  of  (dW/dV)m.n  with  time  for  C 
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